#include "fsw.h"

static FILE *fh_load, *fh_save;
static char runname[80], load_name[80], save_name[80];
static char param[160], value[160], line[160];
static fftw_complex *tmp, *tmp2;

void fload(char** input, ctrl_ptr ptr)
{
  if (input[1] == NULL)
  {
    printf("Usage:\t%s input\n",input[0]);
    exit(1);
  }
  strcpy(runname, input[1]);
  sprintf(load_name,"%s.in",runname);
  sprintf(save_name,"%s.out",runname);
  fh_load = fopen(load_name,"r");
  if (fh_load == NULL) 
  {
    printf("%s\tno such file or file corrupt.\n",load_name);
    exit(1);
  }
  while (fgets(line, 160, fh_load) != NULL)
  {
    sscanf(line, "%s\t%s", param, value);
      if (strcmp(param,"#gridpoints") == 0) ptr->n = atoi(value);
      if (strcmp(param,"#write_skip") == 0) ptr->tskip = atof(value);
      if (strcmp(param,"#box_length") == 0) ptr->l = atof(value);
      if (strcmp(param,"#cutoff_knum") == 0) ptr->kcut = atof(value);
      if (strcmp(param,"#art_viscos") == 0) ptr->av = atof(value);
      if (strcmp(param,"#water_depth")	== 0) ptr->h = atof(value);
      if (strcmp(param,"#grav_accel")	== 0) ptr->g = atof(value);
      if (strcmp(param,"#max_time") == 0)   ptr->t_max = atof(value);
      if (strcmp(param,"#cfl_num") == 0)    ptr->cfl = atof(value);
      if (strcmp(param,"#threads") == 0)    ptr->n_threads = atoi(value);
      if (strcmp(param,"#k_power") == 0)    ptr->k_power = atoi(value);
      if (strcmp(param,"#epsilon") == 0)    ptr->eps = atof(value);
  }
  fclose(fh_load);
  ptr->dx = (ptr->l)/(ptr->n); ptr->dt = (ptr->cfl)*(ptr->dx);
  printf("Read:\n%4d grid points\n%.2f box length\n", ptr->n, ptr->l);
  printf("%.2f water depth\n%.2f gravity acceleration\n",  ptr->h, ptr->g );
  printf("%.2f termination time\n%.2d Wave Number Cut-Off\n", ptr->t_max, ptr->kcut);
  printf("%.2f cfl number\n%d threads\n%d hyper viscosity power\n", ptr->cfl, ptr->n_threads, ptr->k_power);
  ptr->dx = (ptr->l)/(ptr->n);
  ptr->dt = (ptr->cfl)*(ptr->dx);
  ptr->psi = malloc(sizeof(fftw_complex)*(ptr->n));
  ptr->z   = malloc(sizeof(fftw_complex)*(ptr->n));
  ptr->dps = malloc(sizeof(fftw_complex)*(ptr->n));
  ptr->rps = malloc(sizeof(fftw_complex)*(ptr->n));
  ptr->dz  = malloc(sizeof(fftw_complex)*(ptr->n));
  ptr->trpj= malloc(sizeof(fftw_complex)*(ptr->n));
  ptr->fps = malloc(sizeof(fftw_complex)*(ptr->n));
  ptr->fz  = malloc(sizeof(fftw_complex)*(ptr->n));
  fh_load = fopen(load_name,"r"); 
  int i, j = 0;
  char valuex[160],valuepr[160],valuepi[160],valuezr[160],valuezi[160];
  while (fgets(line, 160, fh_load) != NULL)
  {
    sscanf(line,"%s\t\t\t%s\t\t%s\t\t%s\t\t\t%s", valuex, valuepr, valuepi, valuezr, valuezi); 
    j = 0;
    if (strcmp(valuex,"#1.") == 0)
    {
      printf("IC data found. Reading\n");
      char *empty = fgets(line, 160, fh_load);
      while (fgets(line, 160, fh_load) != NULL)
      {
        sscanf(line, "%s\t%s\t%s\t%s\t%s\n", valuex, valuepr, valuepi, valuezr, valuezi);
        if (j < ptr->n) 
	{
          ptr->psi[j] = atof(valuepr) + 1.I*atof(valuepi);
          ptr->z[j]   = atof(valuezr) + 1.I*atof(valuezi);
          j++;
	}
      }
      break;
    } 
  }
  fclose(fh_load);
}

void fsave(char *savename, ctrl_ptr ctrl)
{
  int k;
  fh_save = fopen(savename,"w");
  fprintf(fh_save,"#1.x 2.psi 3.re(z) 4.im(z) 5. spectrum(psi)\n");
  fprintf(fh_save,"#6.spectrum(z) 7.psi_u 8.R psi_u \n"); 
  fprintf(fh_save,"#9.re(dz) 10.im(dz) 11.T{(R psi_u)/J}\n\n");
  /*fprintf(fh_save,"17. re(rhs psi) 18. im(rhs psi) 19. re(rhs z) 18. im(rhs z)\n");*/
  for ( k = 0; k < (ctrl->n); k++ )
  {
    fprintf(fh_save, "%.15e\t%.15e\t%.15e\t%.15e\t%.15e\t", -0.5*(ctrl->l)+(k+0.5)*(ctrl->dx),  creal(ctrl->psi[k]), creal(ctrl->z[k]),  cimag(ctrl->z[k]), cabs(ctrl->fps[k]));
    fprintf(fh_save, "%.15e\t%.15e\t%.15e\t%.15e\t%.15e\t", cabs(ctrl->fz[k]), creal(ctrl->dps[k]), creal(ctrl->rps[k]), creal(ctrl->dz[k]),  cimag(ctrl->dz[k]) ); 
    fprintf(fh_save, "%.15e\n",  creal(ctrl->trpj[k]));
    /*fprintf(fh_save, "%.15e\t%.15e\t%.15e\t%.15e\n", creal(ctrl->rp1[k]),  cimag(ctrl->rp1[k]), creal(ctrl->rz1[k]),  cimag(ctrl->rz1[k]));*/
  }
  fclose(fh_save);
}
